library(car)
library(mosaic)
library(DT)
library(tidyverse)

Read in the Data

files <- dir()
rdat <- read.csv(grep(".csv", files, value=TRUE), header=TRUE)
as_tibble(rdat)
# A tibble: 2,600 × 11
        y    x1    x2    x3    x4    x5    x6    x7      x8    x9   x10
    <dbl> <int> <int> <int> <int> <int> <int> <int>   <dbl> <int> <int>
 1 -4.38      0     0     0     0     0     0     0 -0.0795     1    -1
 2  7.66      1     1     0     0     0     1     0 -2.54       1    -1
 3 -6.82      1     1     0     0     0     1     0 -2.36      -1    -1
 4 -1.23      1     1     0     0     0     1     0  1.45      -1    -1
 5  4.01      0     0     0     0     0     0     0  0.690      1     1
 6  0.730     0     0     0     0     0     0     0  0.175      1    -1
 7 -6.80      1     1     0     0     0     1     0  1.72       1     1
 8  4.99      0     0     0     0     0     0     0 -0.419     -1     1
 9  5.70      1     1     0     0     0     1     0 -0.199      1    -1
10 -1.98      0     0     0     0     0     0     0  0.302     -1    -1
# ℹ 2,590 more rows

Create first Pairs Plot

pairs(rdat, pch=16, cex=1, panel=panel.smooth, col=rgb(.8,.8,.8))

Start with everything except x1 as x1= x2

rdat <- rdat[,-2]
lm1 <- lm(y ~ ., data=rdat)
summary(lm1)
## 
## Call:
## lm(formula = y ~ ., data = rdat)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -34.189  -4.180   0.037   4.136  32.647 
## 
## Coefficients: (4 not defined because of singularities)
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept) -0.23119    0.24247  -0.953  0.34044   
## x2           0.17530    0.35238   0.497  0.61890   
## x3           1.33716    0.46347   2.885  0.00395 **
## x4                NA         NA      NA       NA   
## x5                NA         NA      NA       NA   
## x6                NA         NA      NA       NA   
## x7                NA         NA      NA       NA   
## x8           0.31343    0.13045   2.403  0.01634 * 
## x9           0.07266    0.16075   0.452  0.65131   
## x10          0.14065    0.16086   0.874  0.38201   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.192 on 2594 degrees of freedom
## Multiple R-squared:  0.00613,    Adjusted R-squared:  0.004214 
## F-statistic:   3.2 on 5 and 2594 DF,  p-value: 0.006962

Sample size is large, try for more.

lm2 <- lm(y ~ .^5, data=rdat)
summary(lm2)$coefficients %>%
  as.data.frame() %>%
  arrange(`Pr(>|t|)`) %>%
  datatable()

Pull out significant terms

lm3 <- lm(y ~  x8 +  x3: x8 +  x2: x8: x9: x10 +  x8: x9: x10 +  x3: x8: x9: x10, data=rdat)
summary(lm3)
## 
## Call:
## lm(formula = y ~ x8 + x3:x8 + x2:x8:x9:x10 + x8:x9:x10 + x3:x8:x9:x10, 
##     data = rdat)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -30.8110  -4.1093  -0.0493   4.2569  25.1425 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -0.08016    0.12479  -0.642  0.52071    
## x8            0.31071    0.11487   2.705  0.00688 ** 
## x8:x3        -0.49795    0.24344  -2.045  0.04091 *  
## x8:x9:x10    -6.89563    0.23465 -29.387  < 2e-16 ***
## x8:x2:x9:x10  7.22481    0.26916  26.842  < 2e-16 ***
## x8:x3:x9:x10 13.20240    0.31808  41.507  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.355 on 2594 degrees of freedom
## Multiple R-squared:  0.402,  Adjusted R-squared:  0.4009 
## F-statistic: 348.8 on 5 and 2594 DF,  p-value: < 2.2e-16
palette(c("skyblue","orange","green","purple","steelblue","red","green3","black"))
pairs(cbind(R=lm3$res, Fit=lm3$fit, rdat), pch=16, cex=1, panel=panel.smooth, col=interaction(rdat$x2, rdat$x3, rdat$x9, rdat$x10, drop=TRUE))

Zoom in

ggplot(rdat, aes( x= x8, y=y, color=interaction( x2, x3, x9, x10))) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

ggplot(rdat, aes( x= x8, y=y, color=interaction( x2, x3, x9, x10, x4))) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

ggplot(rdat, aes( x= x8, y=y, color=interaction( x2, x3, x9, x10, x5))) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

ggplot(rdat, aes( x= x8, y=y, color=interaction( x2, x3, x9, x10, x6))) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

ggplot(rdat, aes( x= x8, y=y, color=interaction( x2, x3, x9, x10, x7))) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

ggplot(rdat, aes( x= x8, y=y, color=interaction( x2, x3, x9, x10))) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

Try my own switches

rdat2 <- rdat %>%
  mutate(mycolor = case_when(
     x2==0 &  x3==0 &  x9==-1 &  x10==-1 ~ "Cubic Down",
     x2==0 &  x3==0 &  x9==1 &  x10==1 ~ "Cubic Down",
     x2==0 &  x3==0 &  x9==1 &  x10==-1 ~ "Cubic Up",
     x2==0 &  x3==0 &  x9==-1 &  x10==1 ~ "Cubic Up",
     x2==1 &  x3==0 &  x9==-1 &  x10==-1 ~ "Quintic Up",
     x2==1 &  x3==0 &  x9==1 &  x10==1 ~ "Quintic Up",
     x2==1 &  x3==0 &  x9==1 &  x10==-1 ~ "Quintic Down",
     x2==1 &  x3==0 &  x9==-1 &  x10==1 ~ "Quintic Down",
     x2==0 &  x3==1 &  x9==-1 &  x10==-1  &  x8>0 ~ "Quadratic Up Right",
     x2==0 &  x3==1 &  x9==1 &  x10==1 &  x8>0~ "Quadratic Up Right",
     x2==0 &  x3==1 &  x9==1 &  x10==-1  &  x8<0 ~ "Quadratic Up Left",
     x2==0 &  x3==1 &  x9==-1 &  x10==1 &  x8<0~ "Quadratic Up Left",
     x2==0 &  x3==1 &  x9==-1 &  x10==1 &  x8>0~ "Quadratic Down Right",
     x2==0 &  x3==1 &  x9==1 &  x10==-1 &  x8>0~ "Quadratic Down Right",
     x2==0 &  x3==1 &  x9==-1 &  x10==-1 &  x8<0~ "Quadratic Down Left",
     x2==0 &  x3==1 &  x9==1 &  x10==1 &  x8<0~ "Quadratic Down Left",
  ))

ggplot(rdat2, aes( x= x8, y=y, color=mycolor)) + 
  geom_point() + 
  facet_wrap(~interaction( x2, x3, x9, x10))

ggplot(rdat2, aes( x= x8, y=y, color=mycolor)) + 
  geom_point() + 
  facet_wrap(~mycolor)

lm4 <- lm(y ~ mycolor +  x8:mycolor + I( x8^2):mycolor + I( x8^3):mycolor + I( x8^4):mycolor + I( x8^5):mycolor, data=rdat2)
summary(lm4)
## 
## Call:
## lm(formula = y ~ mycolor + x8:mycolor + I(x8^2):mycolor + I(x8^3):mycolor + 
##     I(x8^4):mycolor + I(x8^5):mycolor, data = rdat2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12.2357  -2.5761  -0.0648   2.6062  13.4635 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           -0.14466    0.32037  -0.452  0.65164    
## mycolorCubic Up                        0.05787    0.44283   0.131  0.89603    
## mycolorQuadratic Down Left           -36.98867   24.74862  -1.495  0.13515    
## mycolorQuadratic Down Right          -70.42413   16.46831  -4.276 1.97e-05 ***
## mycolorQuadratic Up Left             -15.26175   24.37863  -0.626  0.53135    
## mycolorQuadratic Up Right             53.60584   17.47839   3.067  0.00219 ** 
## mycolorQuintic Down                    0.04223    0.45693   0.092  0.92636    
## mycolorQuintic Up                     -0.51824    0.45716  -1.134  0.25706    
## mycolorCubic Down:x8                  11.58827    0.93520  12.391  < 2e-16 ***
## mycolorCubic Up:x8                    -9.30418    0.92190 -10.092  < 2e-16 ***
## mycolorQuadratic Down Left:x8        -28.64204  130.31594  -0.220  0.82605    
## mycolorQuadratic Down Right:x8       228.28264   93.11730   2.452  0.01429 *  
## mycolorQuadratic Up Left:x8         -190.93491  123.56317  -1.545  0.12241    
## mycolorQuadratic Up Right:x8        -126.48696   94.55749  -1.338  0.18112    
## mycolorQuintic Down:x8                -2.68016    0.50609  -5.296 1.29e-07 ***
## mycolorQuintic Up:x8                   3.12033    0.50134   6.224 5.65e-10 ***
## mycolorCubic Down:I(x8^2)             -0.24955    1.06256  -0.235  0.81434    
## mycolorCubic Up:I(x8^2)                0.13598    1.00547   0.135  0.89243    
## mycolorQuadratic Down Left:I(x8^2)    79.11401  248.27532   0.319  0.75001    
## mycolorQuadratic Down Right:I(x8^2) -317.74366  187.34056  -1.696  0.08999 .  
## mycolorQuadratic Up Left:I(x8^2)    -398.50591  228.15814  -1.747  0.08082 .  
## mycolorQuadratic Up Right:I(x8^2)    116.86690  183.32486   0.637  0.52387    
## mycolorQuintic Down:I(x8^2)            0.08022    0.31566   0.254  0.79940    
## mycolorQuintic Up:I(x8^2)              0.42432    0.29088   1.459  0.14477    
## mycolorCubic Down:I(x8^3)            -18.30246    1.97640  -9.261  < 2e-16 ***
## mycolorCubic Up:I(x8^3)               13.94968    1.90453   7.324 3.20e-13 ***
## mycolorQuadratic Down Left:I(x8^3)   128.44133  217.91564   0.589  0.55564    
## mycolorQuadratic Down Right:I(x8^3)  229.56497  171.47661   1.339  0.18077    
## mycolorQuadratic Up Left:I(x8^3)    -334.51958  194.65175  -1.719  0.08582 .  
## mycolorQuadratic Up Right:I(x8^3)    -48.29833  161.95525  -0.298  0.76556    
## mycolorQuintic Down:I(x8^3)            2.88979    0.30231   9.559  < 2e-16 ***
## mycolorQuintic Up:I(x8^3)             -2.88273    0.28754 -10.025  < 2e-16 ***
## mycolorCubic Down:I(x8^4)              0.32585    0.64566   0.505  0.61383    
## mycolorCubic Up:I(x8^4)               -0.08951    0.59065  -0.152  0.87956    
## mycolorQuadratic Down Left:I(x8^4)    65.96118   89.51113   0.737  0.46125    
## mycolorQuadratic Down Right:I(x8^4)  -83.82155   72.65647  -1.154  0.24874    
## mycolorQuadratic Up Left:I(x8^4)    -126.24773   77.80473  -1.623  0.10479    
## mycolorQuadratic Up Right:I(x8^4)      7.53986   66.27789   0.114  0.90944    
## mycolorQuintic Down:I(x8^4)            0.01270    0.05255   0.242  0.80901    
## mycolorQuintic Up:I(x8^4)             -0.05789    0.04616  -1.254  0.20988    
## mycolorCubic Down:I(x8^5)              1.57592    0.93894   1.678  0.09339 .  
## mycolorCubic Up:I(x8^5)                0.34243    0.87828   0.390  0.69665    
## mycolorQuadratic Down Left:I(x8^5)    11.79943   13.92928   0.847  0.39702    
## mycolorQuadratic Down Right:I(x8^5)   11.71172   11.55545   1.014  0.31091    
## mycolorQuadratic Up Left:I(x8^5)     -18.31457   11.78456  -1.554  0.12028    
## mycolorQuadratic Up Right:I(x8^5)      0.35124   10.19107   0.034  0.97251    
## mycolorQuintic Down:I(x8^5)           -0.46370    0.03995 -11.607  < 2e-16 ***
## mycolorQuintic Up:I(x8^5)              0.45323    0.03624  12.507  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.922 on 2552 degrees of freedom
## Multiple R-squared:  0.7759, Adjusted R-squared:  0.7718 
## F-statistic:   188 on 47 and 2552 DF,  p-value: < 2.2e-16
ggplot(rdat2, aes( x= x8, y=y, color=mycolor)) + 
  geom_point(pch=1) + 
  geom_point(aes(y=lm4$fit), cex=0.5) + 
  facet_wrap(~mycolor)

Reduce

Played with this for quite a while. There must be a beautiful combination of x2, x3, x8, x9 that made these all happen with your original switches. But it was easier just to make my own switches.

Well done.

rdat2 <- rdat2 %>%
  mutate( x11 = ifelse(mycolor=="Quadratic Down Right", 1, 0),
          x12 = ifelse(mycolor=="Quadratic Up Right", 1, 0),
          x13 = ifelse(mycolor=="Cubic Down", 1, 0),
          x14 = ifelse(mycolor=="Cubic Up", 1, 0),
          x15 = ifelse(mycolor=="Quintic Down", 1, 0),
          x16 = ifelse(mycolor=="Quintic Up", 1, 0),
          x17 = ifelse(mycolor=="Quadratic Down Left", 1, 0),
          x18 = ifelse(mycolor=="Quadratic Up Left", 1, 0))

lm5 <- lm(y ~  x11 +  x11: x8 +   x11:I( x8^2) + #Quadratic Down Right
             x12 +    x12: x8 +  x12:I( x8^2) +  #Quadratic Up Right
             x17 +  x17: x8 +  x17:I( x8^2) +   #Quadratic Down Left
             x18 +  x18: x8 +  x18:I( x8^2) +   #Quadratic Up Left
             x13: x8 +   x13:I( x8^3) +        #Cubic Down
             x14: x8 +   x14:I( x8^3) +        #Cubic Up
             x15: x8 +   x15:I( x8^3) +    x15:I( x8^5) +     #Quintic Down
             x16: x8 +   x16:I( x8^3) +    x16:I( x8^5)     #Quintic Up
            , data=rdat2)
summary(lm5)
## 
## Call:
## lm(formula = y ~ x11 + x11:x8 + x11:I(x8^2) + x12 + x12:x8 + 
##     x12:I(x8^2) + x17 + x17:x8 + x17:I(x8^2) + x18 + x18:x8 + 
##     x18:I(x8^2) + x13:x8 + x13:I(x8^3) + x14:x8 + x14:I(x8^3) + 
##     x15:x8 + x15:I(x8^3) + x15:I(x8^5) + x16:x8 + x16:I(x8^3) + 
##     x16:I(x8^5), data = rdat2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12.0607  -2.6280  -0.0454   2.6323  13.0411 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -0.05800    0.08474  -0.684    0.494    
## x11         -36.12331    1.72245 -20.972  < 2e-16 ***
## x12          35.41434    1.89730  18.666  < 2e-16 ***
## x17         -34.54107    2.67490 -12.913  < 2e-16 ***
## x18          37.99485    2.35524  16.132  < 2e-16 ***
## x11:x8       56.23008    3.19009  17.626  < 2e-16 ***
## x11:I(x8^2) -22.77164    1.26613 -17.985  < 2e-16 ***
## x8:x12      -52.77655    3.24250 -16.277  < 2e-16 ***
## I(x8^2):x12  21.06233    1.21148  17.386  < 2e-16 ***
## x8:x17      -49.97923    4.33837 -11.520  < 2e-16 ***
## I(x8^2):x17 -19.87236    1.61987 -12.268  < 2e-16 ***
## x8:x18       57.71141    3.81713  15.119  < 2e-16 ***
## I(x8^2):x18  23.23116    1.38806  16.736  < 2e-16 ***
## x8:x13       10.34734    0.52746  19.617  < 2e-16 ***
## x13:I(x8^3) -15.12264    0.43847 -34.490  < 2e-16 ***
## x8:x14       -9.60837    0.51705 -18.583  < 2e-16 ***
## I(x8^3):x14  14.68082    0.41038  35.774  < 2e-16 ***
## x8:x15       -2.64218    0.50620  -5.220 1.94e-07 ***
## I(x8^3):x15   2.85507    0.30207   9.452  < 2e-16 ***
## x15:I(x8^5)  -0.45953    0.03993 -11.509  < 2e-16 ***
## x8:x16        3.19876    0.50060   6.390 1.96e-10 ***
## I(x8^3):x16  -2.92284    0.28742 -10.169  < 2e-16 ***
## I(x8^5):x16   0.45776    0.03626  12.625  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.934 on 2577 degrees of freedom
## Multiple R-squared:  0.7723, Adjusted R-squared:  0.7704 
## F-statistic: 397.3 on 22 and 2577 DF,  p-value: < 2.2e-16
ggplot(rdat2, aes( x= x8, y=y, color=mycolor)) + 
  geom_point(pch=1) + 
  geom_point(aes(y=lm5$fit), cex=0.5) + 
  facet_wrap(~mycolor)

Final Guess

lm

rdat2 <- rdat %>%
  mutate(mycolor = case_when(
     x2==0 &  x3==0 &  x9==-1 &  x10==-1 ~ "Cubic Down",
     x2==0 &  x3==0 &  x9==1 &  x10==1 ~ "Cubic Down",
     x2==0 &  x3==0 &  x9==1 &  x10==-1 ~ "Cubic Up",
     x2==0 &  x3==0 &  x9==-1 &  x10==1 ~ "Cubic Up",
     x2==1 &  x3==0 &  x9==-1 &  x10==-1 ~ "Quintic Up",
     x2==1 &  x3==0 &  x9==1 &  x10==1 ~ "Quintic Up",
     x2==1 &  x3==0 &  x9==1 &  x10==-1 ~ "Quintic Down",
     x2==1 &  x3==0 &  x9==-1 &  x10==1 ~ "Quintic Down",
     x2==0 &  x3==1 &  x9==-1 &  x10==-1  &  x8>0 ~ "Quadratic Up Right",
     x2==0 &  x3==1 &  x9==1 &  x10==1 &  x8>0~ "Quadratic Up Right",
     x2==0 &  x3==1 &  x9==1 &  x10==-1  &  x8<0 ~ "Quadratic Up Left",
     x2==0 &  x3==1 &  x9==-1 &  x10==1 &  x8<0~ "Quadratic Up Left",
     x2==0 &  x3==1 &  x9==-1 &  x10==1 &  x8>0~ "Quadratic Down Right",
     x2==0 &  x3==1 &  x9==1 &  x10==-1 &  x8>0~ "Quadratic Down Right",
     x2==0 &  x3==1 &  x9==-1 &  x10==-1 &  x8<0~ "Quadratic Down Left",
     x2==0 &  x3==1 &  x9==1 &  x10==1 &  x8<0~ "Quadratic Down Left",
  )) %>%
  mutate( x11 = ifelse(mycolor=="Quadratic Down Right", 1, 0),
          x12 = ifelse(mycolor=="Quadratic Up Right", 1, 0),
          x13 = ifelse(mycolor=="Cubic Down", 1, 0),
          x14 = ifelse(mycolor=="Cubic Up", 1, 0),
          x15 = ifelse(mycolor=="Quintic Down", 1, 0),
          x16 = ifelse(mycolor=="Quintic Up", 1, 0),
          x17 = ifelse(mycolor=="Quadratic Down Left", 1, 0),
          x18 = ifelse(mycolor=="Quadratic Up Left", 1, 0))

final.lm <- lm(y ~  x11 +  x11: x8 +   x11:I( x8^2) + #Quadratic Down Right
             x12 +    x12: x8 +  x12:I( x8^2) +  #Quadratic Up Right
             x17 +  x17: x8 +  x17:I( x8^2) +   #Quadratic Down Left
             x18 +  x18: x8 +  x18:I( x8^2) +   #Quadratic Up Left
             x13: x8 +   x13:I( x8^3) +        #Cubic Down
             x14: x8 +   x14:I( x8^3) +        #Cubic Up
             x15: x8 +   x15:I( x8^3) +    x15:I( x8^5) +     #Quintic Down
             x16: x8 +   x16:I( x8^3) +    x16:I( x8^5)     #Quintic Up
            , data=rdat2)
summary(final.lm)

Call:
lm(formula = y ~ x11 + x11:x8 + x11:I(x8^2) + x12 + x12:x8 + 
    x12:I(x8^2) + x17 + x17:x8 + x17:I(x8^2) + x18 + x18:x8 + 
    x18:I(x8^2) + x13:x8 + x13:I(x8^3) + x14:x8 + x14:I(x8^3) + 
    x15:x8 + x15:I(x8^3) + x15:I(x8^5) + x16:x8 + x16:I(x8^3) + 
    x16:I(x8^5), data = rdat2)

Residuals:
     Min       1Q   Median       3Q      Max 
-12.0607  -2.6280  -0.0454   2.6323  13.0411 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -0.05800    0.08474  -0.684    0.494    
x11         -36.12331    1.72245 -20.972  < 2e-16 ***
x12          35.41434    1.89730  18.666  < 2e-16 ***
x17         -34.54107    2.67490 -12.913  < 2e-16 ***
x18          37.99485    2.35524  16.132  < 2e-16 ***
x11:x8       56.23008    3.19009  17.626  < 2e-16 ***
x11:I(x8^2) -22.77164    1.26613 -17.985  < 2e-16 ***
x8:x12      -52.77655    3.24250 -16.277  < 2e-16 ***
I(x8^2):x12  21.06233    1.21148  17.386  < 2e-16 ***
x8:x17      -49.97923    4.33837 -11.520  < 2e-16 ***
I(x8^2):x17 -19.87236    1.61987 -12.268  < 2e-16 ***
x8:x18       57.71141    3.81713  15.119  < 2e-16 ***
I(x8^2):x18  23.23116    1.38806  16.736  < 2e-16 ***
x8:x13       10.34734    0.52746  19.617  < 2e-16 ***
x13:I(x8^3) -15.12264    0.43847 -34.490  < 2e-16 ***
x8:x14       -9.60837    0.51705 -18.583  < 2e-16 ***
I(x8^3):x14  14.68082    0.41038  35.774  < 2e-16 ***
x8:x15       -2.64218    0.50620  -5.220 1.94e-07 ***
I(x8^3):x15   2.85507    0.30207   9.452  < 2e-16 ***
x15:I(x8^5)  -0.45953    0.03993 -11.509  < 2e-16 ***
x8:x16        3.19876    0.50060   6.390 1.96e-10 ***
I(x8^3):x16  -2.92284    0.28742 -10.169  < 2e-16 ***
I(x8^5):x16   0.45776    0.03626  12.625  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3.934 on 2577 degrees of freedom
Multiple R-squared:  0.7723,    Adjusted R-squared:  0.7704 
F-statistic: 397.3 on 22 and 2577 DF,  p-value: < 2.2e-16

Math Model

\[ Y_i = \beta_{0} + \beta_{1} x_{11} + \beta_{2} x_{12} + \beta_{3} x_{17} + \beta_{4} x_{18} + \beta_{5} x_{11} x_{8} + \beta_{6} x_{11} x_{8^2} + \beta_{7} x_{8} x_{12} + \beta_{8} x_{8^2} x_{12} + \beta_{9} x_{8} x_{17} + \beta_{10} x_{8^2} x_{17} + \beta_{11} x_{8} x_{18} + \beta_{12} x_{8^2} x_{18} + \beta_{13} x_{8} x_{13} + \beta_{14} x_{13} x_{8^3} + \beta_{15} x_{8} x_{14} + \beta_{16} x_{8^3} x_{14} + \beta_{17} x_{8} x_{15} + \beta_{18} x_{8^3} x_{15} + \beta_{19} x_{15} x_{8^5} + \beta_{20} x_{8} x_{16} + \beta_{21} x_{8^3} x_{16} + \beta_{22} x_{8^5} x_{16} + \epsilon_i \]

Plot

palette(c("skyblue","orange","green","purple","steelblue","red","green3","black"))

plot(y ~  x8, data=rdat2, col=as.factor(mycolor))
points(final.lm$fit ~  x8, data=rdat2, col=as.factor(mycolor), pch=16, cex=0.5)

b <- coef(final.lm)

drawit <- function( x11=0,  x12=0,  x13=0,  x14=0,  x15=0,  x16=0,  x17=0,  x18=0, i=1){
  curve(b[1] + b[2]* x11 + b[3]* x12 + b[4]* x17 + b[5]* x18 + b[6]* x11* x8 + b[7]* x11* x8^2 + b[8]* x8* x12 + b[9]* x8^2* x12 + b[10]* x8* x17 + b[11]* x8^2* x17 + b[12]* x8* x18 + b[13]* x8^2* x18 + b[14]* x8* x13 + b[15]* x13* x8^3 + b[16]* x8* x14 + b[17]* x8^3* x14 + b[18]* x8* x15 + b[19]* x8^3* x15 + b[20]* x15* x8^5 + b[21]* x8* x16 + b[22]* x8^3* x16 + b[23]* x8^5* x16, add=TRUE,  xname="x8", col=palette()[i])  
}

drawit(1,0,0,0,0,0,0,0,4)
drawit(0,1,0,0,0,0,0,0,6)
drawit(0,0,1,0,0,0,0,0,1)
drawit(0,0,0,1,0,0,0,0,2)
drawit(0,0,0,0,1,0,0,0,7)
drawit(0,0,0,0,0,1,0,0,8)
drawit(0,0,0,0,0,0,1,0,3)
drawit(0,0,0,0,0,0,0,1,5)

with(rdat2, levels(interaction( x11, x12, x13, x14, x15, x16, x17, x18, drop=TRUE)))
## [1] "1.0.0.0.0.0.0.0" "0.1.0.0.0.0.0.0" "0.0.1.0.0.0.0.0" "0.0.0.1.0.0.0.0"
## [5] "0.0.0.0.1.0.0.0" "0.0.0.0.0.1.0.0" "0.0.0.0.0.0.1.0" "0.0.0.0.0.0.0.1"